/* paulslib.h
 *
 * Copyright 2006-2012 David G. Barnes, Paul Bourke, Christopher Fluke
 *
 * This file is part of S2PLOT.
 *
 * S2PLOT is free software: you can redistribute it and/or modify it
 * under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * S2PLOT is distributed in the hope that it will be useful, but
 * WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with S2PLOT.  If not, see <http://www.gnu.org/licenses/>. 
 *
 * We would appreciate it if research outcomes using S2PLOT would
 * provide the following acknowledgement:
 *
 * "Three-dimensional visualisation was conducted with the S2PLOT
 * progamming library"
 *
 * and a reference to
 *
 * D.G.Barnes, C.J.Fluke, P.D.Bourke & O.T.Parry, 2006, Publications
 * of the Astronomical Society of Australia, 23(2), 82-93.
 *
 * $Id: paulslib.h 5786 2012-10-09 01:47:31Z dbarnes $
 *
 */

#ifndef PAULSLIB_H
#define PAULSLIB_H

#if defined(__cplusplus) && !defined(S2_CPPBUILD)
extern "C" {
#endif

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#ifndef WIN32
#include <sys/time.h>
#include <sys/resource.h>
#include <termios.h>
#include <unistd.h>
#else
#include <windows.h>
#include <mmsystem.h>
#endif

#ifndef ABS
#define ABS(x) (x < 0 ? -(x) : (x))
#endif
  
#ifndef MIN
#define MIN(x,y) (x < y ? x : y)
#endif

#ifndef MAX
#define MAX(x,y) (x > y ? x : y)
#endif

#define SIGN(x) (x < 0 ? (-1) : 1)
#define MODULUS(p) (sqrt(p.x*p.x + p.y*p.y + p.z*p.z))
#define DISTANCE(p1,p2) \
	(sqrt((p1.x-p2.x)*(p1.x-p2.x) + \
	(p1.y-p2.y)*(p1.y-p2.y) + \
	(p1.z-p2.z)*(p1.z-p2.z)))
#define CROSSPROD(p1,p2,p3) \
	p3.x = p1.y*p2.z - p1.z*p2.y; \
	p3.y = p1.z*p2.x - p1.x*p2.z; \
	p3.z = p1.x*p2.y - p1.y*p2.x
#define DOTPRODUCT(v1,v2) ( v1.x*v2.x + v1.y*v2.y + v1.z*v2.z )

#define SWAP_2(x) ( (((x) & 0xff) << 8) | ((unsigned short)(x) >> 8) )
#define SWAP_4(x) ( ((x) << 24) | \
         (((x) << 8) & 0x00ff0000) | \
         (((x) >> 8) & 0x0000ff00) | \
         ((x) >> 24) )
#define FIX_SHORT(x) (*(unsigned short *)&(x) = SWAP_2(*(unsigned short *)&(x)))
#define FIX_INT(x)   (*(unsigned int *)&(x)   = SWAP_4(*(unsigned int *)&(x)))
#define FIX_LONG(x)  FIX_INT(x)
#define FIX_FLOAT(x) FIX_INT(x)

  typedef struct {
    double real,imag;
  } pCOMPLEX;

  typedef struct {
    double x,y;
  } XY;
  
#ifndef XYZSTRUCT
  typedef struct {
#if defined(S2TRIPLEFLOAT)
    float x,y,z;
#else
    double x,y,z;
#endif
  } XYZ;
#define XYZSTRUCT 
#endif
  
  typedef struct {
    float x,y,z;
  } fXYZ;
  
  typedef struct {
    XYZ p1,p2,p3;		/* Vertices */
    XYZ c;         	/* Centroid */
    XYZ n;				/* Normal   */
  } FACET3;
  typedef struct {
    XYZ p1,p2,p3;     /* Vertices */
  } TRIFACET;
  typedef struct {
    XYZ p[3];
  } TRIFACE;
  
  typedef struct PLANE {
    double a,b,c,d;
  } PLANE;
  typedef struct {
    XYZ n;
    float D;
  } PLANE2;
  
  typedef struct {
    XYZ p[3];			/* Vertices */
    XYZ c;				/* Centroid */
    XYZ n[3];			/* Normal   */
  } TRIANGLE;
  typedef struct {
    int p1,p2,p3;
  } ITRIANGLE;
  
  typedef struct {
    XYZ p1,p2,p3,p4;  /* Vertices */
    XYZ c;        	 	/* Centroid */
    XYZ n;         	/* Normal   */
  } FACET4;
  
  typedef struct {
    XYZ p1,p2;
  } EDGE;
  typedef struct {
    int p1,p2;
  } IEDGE;
  
  typedef struct {
    XYZ center;
    double radius;
  } SPHERE;
  
  typedef struct {
    int np;
    XYZ p[4];
  } POLY34;
  
  typedef struct {
    XYZ p[8];
    XYZ n[8];
    double val[8];
  } GRIDCELL;
  
#ifndef COLOURSTRUCT
  typedef struct {
#if defined(S2TRIPLEFLOAT)
    float r,g,b,a;
#else
    double r,g,b;
#endif
  } COLOUR;
#define COLOURSTRUCT
#endif
  
  typedef struct {
    float r,g,b; 
  } fCOLOUR;
  typedef struct {
    double h,s,l;
  } HSL;
  typedef struct {
    double h,s,v;
  } HSV;
  typedef struct {
    double y,cr,cb;
  } YCC;
  typedef struct {
    double r,g,b,a;
  } RGBA;
  typedef struct {
    unsigned char r,g,b;
  } PIXEL;
  typedef struct {
    unsigned char r,g,b,a;
  } PIXELA;
  
  typedef struct {
    XYZ p;
    COLOUR c;
  } CXYZ;
  
#ifndef YES
#define YES 1
#endif

#ifndef NO
#define NO  0
#endif

#ifndef TRUE
#define TRUE  1
#define FALSE 0
#endif

#define CLOCKWISE 1
#define ANTICLOCKWISE (-1)
#define COUNTERCLOCKWISE (-1)
#define CONVEX 1
#define CONCAVE (-1)
#define ESC 27

#define THREEPI         9.424777960769379715387930
#define EPSILON         0.000001
#define EPS             0.00001

#define DTOR            0.0174532925
#define RTOD            57.2957795
#define TWOPI           6.283185307179586476925287
#define PI              3.141592653589793238462643
#define PID2            1.570796326794896619231322
#define PID3            1.047197551196397746154214
#define PID4            0.7853981633974483096156608
#define PIPI            9.869604401089358618834491
#define E               2.718281828459045235360287
#define GAMMA           0.5772156649015328606065120
#define LOG2            0.3010299956639811952137389
#define LOG3            0.4772112547196624372950279
#define LOGPI           0.4971498726941338543512683
#define LOGE            0.4342944819032518276511289
#define LOG2E           1.4426950408889634073
#define LN2             0.6931471805599453094172321
#define LN3             1.098612288668109691395245
#define LNPI            1.144729885849400174143427
#define LNSQRT2         0.3465735902799726547
#define LNGAMMA         -0.549539293981644822337662
#define SQRT2           1.4142135623730950488
#define SQRT3           1.7320508075688772935
#define SQRT5           2.2360679774997896964
#define SQRT7           2.6457513110645905905
#define SQRTPI          1.772453850905516027298167

/* Converts from rand() to range 0-1 ie: 1/(2^15-1) */
#define RANDSCALE       0.0000305185095

/* Prototypes  ---------------------------------------------------------*/
void   SplineCurve(XYZ *,int,int *,int,XYZ *,int);
void   SplineKnots(int *,int,int);
void   SplinePoint(int *,int,int,double,XYZ *,XYZ *);
double SplineBlend(int,int,int *,double);
XYZ    Bezier3(XYZ,XYZ,XYZ,double);
XYZ    Bezier4(XYZ,XYZ,XYZ,XYZ,double);
XYZ    Bezier(XYZ *,int,double);
double BezierBlend(int,double,int);
XYZ    CubicBezier(XYZ,XYZ,XYZ,XYZ,double);
XYZ    VertexInterp(double,XYZ,XYZ,double,double);

double Angle2D(double,double,double,double);
double DotProduct(XYZ,XYZ);
double VectorAngle(XYZ,XYZ);
XYZ    CrossProduct(XYZ,XYZ);
double VectorLength(XYZ,XYZ);
void   SetVectorLength(XYZ *,double);
double Modulus(XYZ);
void   Normalise(XYZ *);
void   VectorPerp(XYZ,XYZ *);
int    VectorEqual(XYZ,XYZ);
XYZ    VectorInvert(XYZ);
XYZ    VectorSub(XYZ,XYZ);
XYZ    VectorAdd(XYZ,XYZ);
XYZ    VectorMul(XYZ,double);
void	 MakeVector(XYZ *,double,double,double);
XY     Vector2DRotate(XY,double);
XYZ    RotateX(XYZ,double);
XYZ    RotateY(XYZ,double);
XYZ    RotateZ(XYZ,double);
XYZ    MidPoint(XYZ,XYZ);
XYZ    Centroid(XYZ *,int);
double FowlerAngle(double,double);
XYZ    CalcNormal(XYZ,XYZ,XYZ);
int    LineLineIntersect(XYZ,XYZ,XYZ,XYZ,XYZ *,XYZ *,double *, double *);
int    FindClosest(XYZ,XYZ *,int,double *,XYZ *);
int    ClipFacet(XYZ *,XYZ,XYZ);
int    ClipLine2Box(XYZ *,XYZ *,XYZ,XYZ);
int    PointInBox(XYZ,XYZ,XYZ);
int    ContourFacet(XYZ *,XYZ,XYZ,XYZ *,XYZ *);
int    LineFacet(XYZ,XYZ,XYZ,XYZ,XYZ,XYZ *);
int    LineIntersect(double,double,double,double,
                     double,double,double,double,
                     double *,double *);
double PointLine2D(XY,XY,XY,XY *,double *);
double PointLine3D(XYZ,XYZ,XYZ,XYZ *,double *);
int    ClockWise(XY *,int);
int    ConvexPolygon(XY *,int);
void   XYZ2Polar(double,double,double,double *,double *,double *);
void   InverseSpherical(XYZ,double *,double *);
void   Polar2XYZ(double,double,double,double,double,double);
int    Triangulate(int,XYZ *,ITRIANGLE *,int *);
int    CircumCircle(double,double,
	       double,double,double,double,double,double,
			 double *,double *,double *);
XYZ    ArbitraryRotate(XYZ,double,XYZ);
XYZ    ArbitraryRotateUnit(XYZ,double,XYZ);
XYZ    ArbitraryRotate2(XYZ,double,XYZ,XYZ);
XYZ    ArbitraryRotate3(XYZ,double,XYZ,XYZ);
double LinearInterpolate(double,double,double);
double CosineInterpolate(double,double,double);
double CubicInterpolate(double,double,double,double,double);
double HermiteInterpolate(double,double,double,double,double,double,double);
void   Contour(double **,int,int,int,int,double *,double *,int,double *,
       void (*drawline)(double,double,double,double,double));
int    RaySphere(XYZ,XYZ,XYZ,double,double *,double *);
int    LinePlane(XYZ,XYZ,PLANE,double *,XYZ *);

/* FILE IO ---------------------------------------------------------*/
int    RLECompress(unsigned char *,unsigned char *,int);
void   RLEUncompress(unsigned char *,unsigned char *,int);
int    ReadLine(FILE *,char *,int);
int    ReadUntil(FILE *,int,char *);
void   ReadToString(FILE *,char *);
void   RightTrim(char *);
int    ReadString(FILE *,char *,int);
int    ReadShort(FILE *,short *,int);
int    WriteShort(FILE *,short,int);
int    ReadUShort(FILE *,unsigned short *,int);
int    WriteUShort(FILE *,unsigned short,int);
int    ReadLong(FILE *,long *,int);
int 	 WriteLong(FILE *,long,int);
int    WriteULong(FILE *,unsigned long,int);
int    ReadInt(FILE *,int *,int);
int    ReadUInt(FILE *,unsigned int *,int);
int    WriteInt(FILE *,int,int);
int    WriteUInt(FILE *,unsigned int,int);
int    ReadDouble(FILE *,double *,int);
int    WriteDouble(FILE *,double,int);
int    ReadFloat(FILE *,float *,int);
int    WriteFloat(FILE *,float,int);
int    ReadSkip(FILE *,int);
void   SkipUntil(FILE *,int);
void   StrToLower(char *);
void   StrToUpper(char *);
void   EchoFile(char *,FILE *);
int    IndexCopy(char *,char *,unsigned int,unsigned int);

/* COLOUR STUFF -----------------------------------------------------*/
COLOUR RampColour2(double,COLOUR,COLOUR);
COLOUR GetColour(double,double,double,int);
int    AsciiColour(double,double,double);
double AbsColour(COLOUR);
void   ClipColour(COLOUR *);
void   ClipRGBA(RGBA *);
void   InvertColour(COLOUR *);
int    EqualColour(COLOUR,COLOUR);
HSL    RGB2HSL(COLOUR);
COLOUR HSL2RGB(HSL);
HSV    RGB2HSV(COLOUR);
COLOUR HSV2RGB(HSV);
YCC    RGB2YCC(COLOUR);
COLOUR YCC2RGB(YCC);

/* POLYGON STUFF -----------------------------------------------------*/
#if defined(S2TRIPLEFLOAT)
void   PlaneParam(XYZ *,float *,float *,float *,double *);
#else
void   PlaneParam(XYZ *, double *, double *, double*, double *);
#endif
XY     PointInPolygon(XY *,double,double);
int    PointIn3DPolygon(int,XYZ *,XYZ,double *);
int    WhichSide(double,double,double,double,double,double);
double PlanePoint(double,double,XYZ,XYZ,XYZ);
int    InTriangle(double,double,XYZ,XYZ,XYZ);
double PolygonLength(XYZ *,int);
double PolygonCurve(XYZ *,int);
int    PolygoniseCube(GRIDCELL,double,TRIANGLE *);
int    PolygoniseTri(GRIDCELL,double,TRIANGLE *,int,int,int,int);
int    PolygonCompare(POLY34,POLY34);
int    CreateNSphere(TRIFACE *,int);

/* GEOMETRIC OUTPUT --------------------------------------------------*/
void   WriteRadFacet(FILE *,XYZ *,int,char *);
void   WriteDXFHeader(FILE *,char *);
void   WriteDXFFooter(FILE *);
void   WriteDXFFacet(FILE *,XYZ *,int);

/* SYSTEM SPECIFIC  --------------------------------------------------*/
double GetRunTime(void);
#ifndef WIN32
void SetKeyboardIO(int);
#endif

/* MISC FUNCTIONS  ---------------------------------------------------*/
int    ClipDouble(double *,double,double);
int    ClipInteger(int *,int,int);
void   SwapDouble(double *,double *);
void   SwapXYZ(XYZ *,XYZ *);
void   MinMaxXYZ(XYZ,XYZ *,XYZ *);
void   StringClean(char *);

double Logarithm(double,double);
long   PowerInt(long,long);
long   IntSqrt(long);
double Factorial(int);
double BesselI0(double);
double Bessel(double,int);
double BesselZero(int,int);

double Determinant22(double a[2][2]);
double Det2x2(double,double,double,double);
int    Inverse22(double a[2][2],double b[2][2]); 

double Determinant33(double a[3][3]);
double Det3x3(double,double,double,double,double,double,double,double,double);
int    Inverse33(double a[3][3],double b[3][3]);

double Determinate44(double a[4][4]);
void   Adjoint44(double a[4][4],double b[4][4]);
int    Inverse44(double a[4][4],double b[4][4]);
void   Identity44(double a[4][4]);
void   Multiply44(double a[4][4],double b[4][4],double c[4][4]);

double Determinant(double **,int);
void CoFactor(double **,int,double **);
void Transpose(double **,int);
void MatrixMul(double **,double **,double **,int);
void WriteMatrix(FILE *,int,double **,int);

/* QUATERNIONS  ---------------------------------------------------*/
typedef struct {
   double r,a,b,c;
} QUATERNION;
QUATERNION Qadd(QUATERNION,QUATERNION);
QUATERNION Qsub(QUATERNION,QUATERNION);
QUATERNION Qmult(QUATERNION,QUATERNION);
QUATERNION Qdiv(QUATERNION,QUATERNION);
QUATERNION Qconj(QUATERNION);
QUATERNION Qinv(QUATERNION);
double Qmod(QUATERNION);

/* SIGNALS  ---------------------------------------------------*/
int    FFT(int,int,double *,double *);
int    FFT2D(pCOMPLEX **, int,int,int);
int    DFT(int,int,double *,double *);
int    Powerof2(int,int *,int *);
void   Correlate(double *,double *,int,double *);
int    GSolve(double **,int,double *);
double AutoCorr(double *,int,int);
double CrossCorr(double *,double *,int,int);
void   Write_AIFF(FILE *,double *,long,int);
void   Write_AU(FILE *,double *,long,int);
void   Write_WAVE(FILE *,float *,int,int);
int    LinRegress(double *,double *,int,double *,double *,double *);

/* pCOMPLEX NUMBERS  -------------------------------------------*/
pCOMPLEX Csqrt(pCOMPLEX);
pCOMPLEX Clog(pCOMPLEX);
pCOMPLEX Cexp(pCOMPLEX);
pCOMPLEX Csin(pCOMPLEX);
pCOMPLEX Ccos(pCOMPLEX);
pCOMPLEX Ctan(pCOMPLEX);
pCOMPLEX Casin(pCOMPLEX);
pCOMPLEX Cacos(pCOMPLEX);
pCOMPLEX Catan(pCOMPLEX);
pCOMPLEX Csinh(pCOMPLEX);
pCOMPLEX Ccosh(pCOMPLEX);
pCOMPLEX Ctanh(pCOMPLEX);
pCOMPLEX Casinh(pCOMPLEX);
pCOMPLEX Cacosh(pCOMPLEX);
pCOMPLEX Catanh(pCOMPLEX);
pCOMPLEX Cadd(pCOMPLEX,pCOMPLEX);
pCOMPLEX Csub(pCOMPLEX,pCOMPLEX);
pCOMPLEX Cmult(pCOMPLEX,pCOMPLEX);
pCOMPLEX Cmultd(pCOMPLEX,double);
pCOMPLEX Cdiv(pCOMPLEX,pCOMPLEX);
pCOMPLEX Cinv(pCOMPLEX);
pCOMPLEX Cdivd(pCOMPLEX,double);
pCOMPLEX Cpowd(pCOMPLEX,double);
double  Cabs(pCOMPLEX);

#if defined(__cplusplus) && !defined(S2_CPPBUILD)
} /* extern "C" { */
#endif

#endif /* PAULSLIB_H */

